#!/usr/bin/env python3.9
# This generates 'analogue' SSTV from a (preferable very high resolution)
# image. It's slow, and Python uses up a ridiculous amount of memory, but
# it's just accessible code for playing with.
# --
# Licenced according to the GNU Public Licence version 3.0
# Chris Baird,, May-June 2022 <vk2cjb@gmail.com>

import math, struct, sys
from PIL import Image

# scanlen is number of samples it needs = 48000 * 195.854ms, etc.
# The "Daytona SSTV Paper" scanline timings seem a bit off... QSSTV is happier
# with these (and probably other decoders using audio sampling at 48k).
# fudge is the 'true' duration of a scanline in audio samples, at least, as
# QSSTV sees it. It's a pain. QSSTV is probably faulty.
# And is QSSTV's handling of PD240 broken?

# scanlen has been adjusted for the optimal display with QSSTV (and hopefully
# other computer-based decoders). The original slowscan specifications being
# in milliseconds is not that compatible with digitizing soundcards and most
# needed ~6 samples added to reduce chromatic aberration. Again, QSSTV isn't
# perfect.

specs = [
    #    name type    VIS scanlen     fudge    modesize
    ( "PD50", 'P',     93,  4394, 18632.85, (320,256)),
    ( "PD90", 'P',     99,  8173, 33748.00, (320,256)),
    ("PD120", 'P',     95,  5838, 24408.65, (640,496)),
    ("PD160", 'P',     98,  9389, 38614.60, (512,400)),
    ("PD180", 'P',     96,  8787, 36206.05, (640,496)),
    ("PD240", 'P',     97, 11736, 48003.25, (640,496)),# decodes poorly in QSSTV
    ("PD290", 'P',     94, 10984, 44992.50, (800,616)),
    ("MP115", 'p', 0x2923, 10704, 43299.00, (320,256)),
    ("MP140", 'p', 0x2a23, 12960, 52323.00, (320,256)),
    ("MP175", 'p', 0x2c23, 16320, 65765.00, (320,256)),
    ( "MR73", 'R', 0x4523,  3320, 27486.68, (320,256)),
    ( "MR90", 'R', 0x4623,  4110, 33823.50, (320,256)),
    ("MR115", 'R', 0x4923,  5288, 43232.20, (320,256)),
    ("MR140", 'R', 0x4a23,  6464, 52640.30, (320,256)),
    ("MR175", 'R', 0x4c23,  8095, 65697.00, (320,256)),
    ("ML180", 'R', 0x8523,  4243, 34879.10, (640,496)),
    ("ML240", 'R', 0x8623,  5683, 46399.80, (640,496)),
    ("ML280", 'R', 0x8923,  6666, 54272.05, (640,496)),
    ("ML320", 'R', 0x8a23,  7626, 61953.15, (640,496))]

idsize = (42,9) # guess wildly what you have to do here...
id = (
    0, 0,0,0,0,0, 0,0, 0,0,0,0,0, 0,0, 0,0,0,0,0, 0,0, 0,0,0,0,0, 0,0, 0,0,0,0,0, 0,0, 0,0,0,0,0, 0,
    0, 1,0,0,0,1, 0,0, 1,0,0,0,1, 0,0, 0,1,1,1,0, 0,0, 0,1,1,1,0, 0,0, 0,0,1,1,1, 0,0, 1,1,1,1,0, 0,
    0, 1,0,0,0,1, 0,0, 1,0,0,1,0, 0,0, 1,0,0,0,1, 0,0, 1,0,0,0,1, 0,0, 0,0,0,0,1, 0,0, 1,0,0,0,1, 0,
    0, 1,0,0,0,1, 0,0, 1,0,1,0,0, 0,0, 0,0,0,0,1, 0,0, 1,0,0,0,0, 0,0, 0,0,0,0,1, 0,0, 1,0,0,0,1, 0,
    0, 1,0,0,0,1, 0,0, 1,1,0,0,0, 0,0, 0,1,1,1,0, 0,0, 1,0,0,0,0, 0,0, 0,0,0,0,1, 0,0, 1,1,1,1,0, 0,
    0, 1,0,0,0,1, 0,0, 1,0,1,0,0, 0,0, 1,0,0,0,0, 0,0, 1,0,0,0,0, 0,0, 0,0,0,0,1, 0,0, 1,0,0,0,1, 0,
    0, 0,1,0,1,0, 0,0, 1,0,0,1,0, 0,0, 1,0,0,0,0, 0,0, 1,0,0,0,1, 0,0, 1,0,0,0,1, 0,0, 1,0,0,0,1, 0,
    0, 0,0,1,0,0, 0,0, 1,0,0,0,1, 0,0, 1,1,1,1,1, 0,0, 0,1,1,1,0, 0,0, 0,1,1,1,0, 0,0, 1,1,1,1,0, 0,
    0, 0,0,0,0,0, 0,0, 0,0,0,0,0, 0,0, 0,0,0,0,0, 0,0, 0,0,0,0,0, 0,0, 0,0,0,0,0, 0,0, 0,0,0,0,0, 0)

altidsize = (36,5)
altid=(0,0,0,0,1,0,1,0,0,1,0,1,0,0,1,1,1,0,0,1,1,1,0,0,0,0,1,0,0,1,1,1,0,0,0,0,
      0,0,0,0,1,0,1,0,0,1,0,1,0,0,0,0,1,0,0,1,0,0,0,0,0,0,1,0,0,1,0,1,0,0,0,0,
      0,0,0,0,1,0,1,0,0,1,1,0,0,0,1,1,1,0,0,1,0,0,0,0,0,0,1,0,0,1,1,1,0,0,0,0,
      0,0,0,0,1,0,1,0,0,1,0,1,0,0,1,0,0,0,0,1,0,0,0,0,1,0,1,0,0,1,0,1,0,0,0,0,
      0,0,0,0,0,1,0,0,0,1,0,1,0,0,1,1,1,0,0,1,1,1,0,0,1,1,1,0,0,1,1,1,0,0,0,0)

modelist = [x[0] for x in specs]
sample_rate, ph, nframes, mode, outfile = 48000, 0.0, 0, -1, ""
greymode, showident, leftident, redident = False, False, False, False
audio, domandel, infile, lastfreq, sine = bytearray(), False, "", 0, 0

###

while len(sys.argv) > 1:
    sys.argv.pop(0)
    s = sys.argv[0]
    if s == "-g": greymode = True
    elif s == "-i": showident = True
    elif s == "-r": redident = True
    elif s == "-l": leftident = True
    elif s == "-m":
        domandel = True
        greymode = True
    elif s == "-o":
        outfile = sys.argv[1]
        sys.argv.pop(0)
    elif s in modelist: mode = modelist.index(s)
    else:
        infile = s
        if not len(infile):
            print ("You needed instructions? oh.")
            sys.exit(1)

if mode == -1:
    print ("No valid mode given.")
    sys.exit(1)

name, modetype, VIS, scanlen, fudge, modesize = specs[mode]
height, identstart = modesize[1], modesize[1] * (11/12)

#

if len(infile):
    im_orig = Image.open(infile)
    width, oheight = im_orig.width, im_orig.height
    if im_orig.mode == "I": depth = 65536.0
    else: depth = 256.0

    ycbcr = im_orig.convert('YCbCr')
    y, cb, cr = ycbcr.split()

    imY = y.resize((width,height)).convert("F")
    imCb = cb.resize((width,height)).convert("F")
    imCr = cr.resize((width,height)).convert("F")
    imCbH = cb.resize((width,height//2)).convert("F")
    imCrH = cr.resize((width,height//2)).convert("F")
else:
    width, oheight, depth = modesize[0] * 32, modesize[1], 65536.0

wantedwidth = (modesize[0] * oheight) / modesize[1]
scanlinescale = (modesize[0] * oheight) / (modesize[1] * scanlen)
leftedge = (width - wantedwidth) / 2

###

def sampleout(freq):
    global ph,nframes,audio,lastfreq,sine
    sine = int(math.sin(ph / sample_rate) * 3276.8) # constant can vary
    if sine < 0: sine += 65536
    audio += bytes([sine & 255,sine >> 8])
    ph += math.tau * freq
    nframes += 1
    lastfreq = freq


def tonegen(freq,samples):
    for _ in range(samples): sampleout(freq)

###

def vis(c):
    tonegen(1900,14400)  # leader - 1900 Hz, 300 ms/570 cycles/14400
    tonegen(1200,480)    # break - 1200 Hz,  10 ms/12 cycles/480
    tonegen(1900,14400)  # leader again
    tonegen(1200,1440)   # startbit - 1200 Hz, 30 ms/36 cycles

    p = 0
    for i in (1,2,4,8,16,32,64,128):
        if i == 128 and (p & 1): c |= 128 # even parity
        if (c & i):
            p += 1
            tonegen(1100,1440)
        else:
            tonegen(1300,1440)

    if c > 255:                 # it's an MMSSTV 16-bit VIS code..
        p = 0
        for i in (256,512,1024,2048,4096,8192,16384,32768):
            if i == 32768 and not (p & 1): c |= 32768 # odd parity
            if (c & i):
                p += 1
                tonegen(1100,1440)
            else:
                tonegen(1300,1440)

    tonegen(1200,1440)          # stop bit - 30 ms

###

def doidentline(y,scanlen,side):
    yt = idsize[0] * (int((12*idsize[1]*y / modesize[1])) % idsize[1])
    for x in range(scanlen):
        xx = int((idsize[0]*x*4 / scanlen)) % idsize[0]
        i = xx + yt
        freq = 1900
        if i < len(id) and id[i]: freq = 2300
        if side and x > (scanlen/4): freq = 1900
        if not side and x < (3*scanlen/4): freq = 1900
        sampleout(freq)

###

def mandel(a,b,n):
    x, y = 0,0
    for i in range(n):
        xx = x*x - y*y + a;
        y = 2*x*y + b;
        x = xx;
        if x*x + y*y > 4: break
    return i


def domandelline(y,scanlen):
    xscale = modesize[0] / (modesize[1] * scanlen)
    for x in range(scanlen):
        # the classic mandelbrot pose
        #x1 = 0.8 - (2.5 * x * xscale)
        #y1 = -1.25 + (2.5 * y / modesize[1])
        #d = 2300 - 145.0*math.log(1+mandel(x1, y1, 300))

        # my favourite..
        x1 = -0.55 + (0.000011 * x * xscale)
        y1 = -0.55 + (0.000011 * y / modesize[1])

        m = mandel(x1, y1, 1500)
        #d = 2300 - 150.0*math.log(1+m)
        #d = 2300 - 800 * (m & 1)
        d = 2300 - 114.3 * (m & 7)

        if d < 1500: d = 1500
        sampleout(d)

###

def getvalue(image,sx,y):
    x = leftedge + (sx * scanlinescale)
    if x < 0 or x >= width: r = 0
    else: r = image.getpixel((x,y)) / depth
    return r


def doline(image,y,scanlen):
    for x in range(scanlen):
        c = getvalue(image,x,y)
        freq = 1500 + (800.0 * c)
        sampleout(freq)

def dodoubleline(image,y,scanlen):
    for x in range(scanlen*2):
        c = getvalue(image,x/2,y)
        freq = 1500 + (800.0 * c)
        sampleout(freq)

###

def pdmode(y,f):
    tonegen(1200,960)
    tonegen(1500,96+f)

    if domandel: domandelline(y, scanlen)
    else: doline(imY, y, scanlen)

    if not greymode:
        doline(imCrH, y/2, scanlen)
        doline(imCbH, y/2, scanlen)

    if greymode:
        if showident and y > identstart:
            if redident:
                doidentline(y,scanlen,leftident)
                tonegen(1900,scanlen)
            else:
                tonegen(1900,scanlen)
                doidentline(y,scanlen,leftident)
        else:
            tonegen(1900,scanlen)
            tonegen(1900,scanlen)

    if domandel: domandelline(y+1, scanlen)
    else: doline(imY, y+1, scanlen)

#

def mrmode(y,f):
    g = [432,5,5,5]
    g[0] += f
    print (g)

    tonegen(1200,g[0])          # (432) sync 9.0 ms
    tonegen(1500,48)            # 1.0ms

    # too much mucking about to include mandelbrot, sorry..
    dodoubleline(imY,y,scanlen)
    tonegen(lastfreq,g[1])

    if not greymode:
        doline(imCr,y,scanlen)
        tonegen(lastfreq,g[2])
        doline(imCb,y,scanlen)
        tonegen(lastfreq,g[3])

    if greymode:
        if showident and y > identstart:
            if redident:
                doidentline(y,scanlen,leftident)
                tonegen(lastfreq,g[2])
                tonegen(1900+g[3],scanlen)
            else:
                tonegen(1900+g[2],scanlen)
                doidentline(y,scanlen,leftident)
                tonegen(lastfreq,g[3])
        else:
            tonegen(1900+g[2],scanlen)
            tonegen(1900+g[3],scanlen)

#

def mpmode(y,f):
    tonegen(1200,432)
    tonegen(1500,48+f)

    if domandel: domandelline(y,scanlen)
    else: doline(imY,y,scanlen)

    if not greymode:
        doline(imCrH,y//2,scanlen)
        doline(imCbH,y//2,scanlen)

    if greymode:
        if showident and y > identstart:
            if redident:
                doidentline(y,scanlen,leftident)
                tonegen(1900,scanlen)
            else:
                tonegen(1900,scanlen)
                doidentline(y,scanlen,leftident)
        else:
            tonegen(1900,scanlen)
            tonegen(1900,scanlen)

    if domandel: domandelline(y+1,scanlen)
    else: doline(imY,y+1,scanlen)

# Main

tonegen(0,24000)
vis(VIS)

if modetype == 'P' or modetype == 'p': ystep = 2
else: ystep = 1
t = 0

for y in range(0,height,ystep):
    print (y, end='\r')
    f = int(y*fudge*0.5 - t)
    onf = nframes

    if modetype == 'P': pdmode(y,f)
    if modetype == 'R': mrmode(y,f)
    if modetype == 'p': mpmode(y,f)

    print ("y=",y, "expectt=",y*fudge*0.5, "realt=", t, nframes - onf, f)
    t += nframes - onf

tonegen(0,24000)

###

if not outfile:
    if domandel: infile = "mandelbrot"
    outfile = infile + "-" + name + ".wav"

wav = open(outfile,"wb")

s = struct.Struct(b"<ccccLccccccccLHHLLHHccccL")
wav.write(s.pack(b'R',b'I',b'F',b'F', # riff_header
                 2*len(audio)+44,     # wav_size
                 b'W',b'A',b'V',b'E', # wave_header
                 b'f',b'm',b't',b' ', # fmt_header
                 16,                  # fmt_chunk_size
                 1,                   # audio_format = PCM
                 1,                   # num_channels
                 sample_rate,
                 sample_rate*1*2,     # byte_rate
                 2,                   # sample_alignment
                 16,                  # bit_depth
                 b'd',b'a',b't',b'a',
                 2*len(audio)))       # data_bytes
wav.write(audio)                      # thanks to bytearray
wav.close()

###
